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1. Introduction 

Scatterometry is a common technique for the characterization of nano-structured surfaces. 
It is an indirect measuring method inferring the properties of the scattering object from the 
light diffracted. The interaction between the scattering object and the incoming electromag- 
netic radiation is simulated numerically. Applying an optimization algorithm, the measured 
diffraction intensities are fitted to the calculated ones to derive the sought properties of 
the sample. It should be noted that the same numerical models are also used to predict, 
e.g. the imaging properties of mask structures in a photolithographic process. A better un- 
derstanding of the accuracy of the optical models is therefore crucial also for many other 
applications than scatterometry. Scatterometry, however, is the natural environment for the 
investigation of the accuracy of optical models because it is explicitly based on the corre- 
lation of geometrical structure and optical properties. Currently, scatterometry is used as 
a relative metrology method for process control and process development. Many attempts 
have been made to establish scatterometry as a traceable and absolute metrological method 
for dimensional measurements of nanostructured surfaces [IHU , some of them including the 
evaluation of measurements on test structures for line roughness [51IE]- 

At the Physikalisch-Technische Bundesanstalt (PTB), detailed scatterometric investiga- 
tions of an EUV photomask with periodic absorber line grating test patterns have been 
performed [318]. It was feasible to derive the information on the line profile by means of rig- 
orous numerical modelling jHEHU] . A detailed analysis of the uncertainty contributions in 
the structure reconstruction process [I] revealed that structure disturbances such as line edge 
or line width roughness have significant impact on the angular distribution of the diffraction 
intensities [H1[T2]. Uncertainties in the actual intensity measurements are less important. 
In [13], a Debye- Waller-type attenuation factor was found using a stochastic model of line 
roughness. It decreases the expectation value of the diffraction intensities perpendicular to 
the absorber lines. The calculations are based upon an analytical model describing the far- 
field diffraction of a disturbed ID binary grating. The Debye- Waller factor found there was 
applied to the scattering of a patterned EUV multilayer photomask, and the roughness values 
found showed a high correlation compared to the rms roughness measured by scanning elec- 
tron microscopy [12]. Within a Monte Carlo simulation of the randomly disturbed EUV line 
grating over several periods, the exponential behaviour of the efficiency attenuation for the 
in-plane scattering has recently been confirmed [H]. However, this recent treatment was for a 
ID structure, which inherently only describes in-plane scatter between the diffraction orders. 
In realistic investigations of roughness in scatterometry [5ll6|[T5] also out-of-plane scatter has 
to be included. It will be shown that in the Fraunhofer far-field approximation, the impact 
of the line roughness is given by a multiplicative term also for 2D structures. While scat- 
terometry basically is the observation of the far-field, the condition for the Fraunhofer and 
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Kirchhoff approximations of the structures being substantially larger than the wavelength is 
not generally met, particularly not for optical scatterometry at state-of-the-art semiconduc- 
tor structures. This condition is equivalent to the existence of several orders of diffraction. 
For the example presented here, EUV scatterometry at structures of about 100 nm critical 
dimensions, it is met. Even for structures smaller than the wavelength, the spatial wave- 
lengths of the perturbations (e.g. stitching errors) may be larger than the wavelength and 
roughness induced scatter as described by the formalism presented here is observable. For 
small structures and small spatial wavelength perturbations, however, no other diffraction 
orders than the zeroth order exist and other descriptions like the effective layer model [T5|[T6] 
must be used to account for roughness. 

For periodic edge roughness, the out-of-plane diffraction is described by Bessel functions. 
This approach can be generalized to random roughness by using a Fourier expansion of the 
rough contours. For the in-plane scatter, the correction factor converges for small rough- 
ness amplitudes to an exponential factor, as derived before in the ID case. The analytical 
approach presented here is also compared with rigorous numerical calculations. The EUV 
photomask used as an example consists of a quarter inch thick substrate with a reflective 
Mo/Si multilayer coating. The multilayer coating is terminated with a protective Si/SiC>2 
capping layer. On top of the multilayer, an absorber stack is deposited, which is then struc- 
tured by e-beam lithography. The paper is structured as follows: in section [2] an analytical 
study of the 2D line roughness is presented. In order to evaluate the scatterometry measure- 
ments we use a FEM based Maxwell solver as simulation tool. In section [3] this numerical 
approach is applied to analyze the scattering caused by line roughness. In section H] the 
results of the two previous sections are compared and discussed. 

2. Analytical modeling 

Within the classical Kirchhoff's diffraction theory, the electromagnetic field is treated as a 
scalar quantity. The grating is considered as an aperture with periodic slits which is irradiated 
by monochromatic radiation. Since the example we are presenting here is a grating consisting 
of periodic absorber lines on a mirror surface, we do not use transmission but reflection. Thus 
the grating is represented by a periodically varied reflectivity r on the sample's surface in 
the (x, y) plane. 
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r(x, y) = ^ 5 (x — jd) * rect ^— ^ for all x G 

j=-n 

Here, "*" stands for the convolution over x. The lines are parallel to the y-axis, cf. Fig. 
[U (a), d is the period of the structure in re-direction, c is the width of the reflecting area. 
For simplicity, we set the reflectivity value of the open reflecting multilayer to 1 and that 
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of the absorber lines to zero, and it is assumed that the length of the grating is infinite 
in y-direction. n G N, and N = In + 1 gives the total number of periods. Henceforth, we 
will refer to this kind of grating, which is invariant in y, as an undisturbed grating and the 
subscript r denotes the quantities referring to the undisturbed grating. In the far field, the 
amplitude E of the electric field is proportional to the Fourier transform of the reflection 
function r given by Eq. (pQ), 

E (k x , k y ) oc 7 {r} (k x , k y ) , (2) 

where J 7 {r} stands for the Fourier transform of r with the wave vectors (k x , k y ) as argu- 
ments. In the scalar approximation, the intensity I r is given by the square of the electrical 
field: 

(gjjl k x c gj^-Q Nk x d \ 
C ^J~ • Id ^5{k y ) . (3) 
2 Sm 2 / 

Eq. ([3]) describes the well-known Fraunhofer diffraction of a grating, whose discrete diffrac- 
tion orders are found at k x = 2nm/d for all m G Z and k y = 0. The corresponding intensity 
values will be called diffraction intensities of order m in x, respectively, throughout the paper. 
5 is the Dirac delta function. We are aware of the simplicity of this binary grating model, 
neglecting the mask's 3D geometry, all material properties, and the reflection characteristics 
of the multilayer coating. Therefore, the form factor given by the sine term in Eq. ([3]) will 
not be found exactly in measurements. The advantage is, however, that based on this model, 
the structure roughness can be treated analytically. In [12j[T3] the impact of the stochastic 
ID roughness was successfully described. The roughness could explain the observed differ- 
ences in the side wall angle as reconstructed by angle-resolved scatterometry and measured 
independently by atomic force microscopy, whereas other possible uncertainty contributions 
could not explain this difference adequately [I] . 

The roughness models we consider are depicted in Fig. [1] (b) and (c). The line edges have 
sinusoidal shapes that are in phase for all edges in the case of line edge roughness (LER), 
keeping a fixed value for the line width. In case of the line width roughness (LWR), the line 
edges are phase shifted by 7r which means that the line centre position is fixed along the lines 
while the line width shows a periodic modulation. The framed boxes in Fig. [1] indicate the unit 
cells calculated with the finite element method (FEM) in section [3] using periodic boundary 
conditions. Note that any random roughness can be expressed by a Fourier expansion as the 
sum of sinusoidal waves. It is therefore sufficient to find a solution for the sinusoidal case and 
the general solution is just a linear superposition. Therefore, this approach is widely used. 
E.g. [T7] presents the numerical modelling of 2D line edge roughness for optical scatterometry 
using the field- stitching method in rigorous coupled-wave analysis. In this paper, we will first 
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derive an analytical description of the line edge modulation and will then compare the results 
to 3D FEM simulations. 

2. A. Line edge roughness 

The reflectivity of a binary line-and-space grating with LER as outlined in Fig. [TJ (b) is given 
by 



X 

a cos (2iry/d r ))) * rect ( — 

c 



for all x G 
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stands for the convolution over x and a cos (2iry/d r ) desribes the sinosoidal 
edge with spatial period d r and amplitude a. In the following equations / stands for the 
LER-disturbed grating instead of r for the undisturbed solution in (J3J). Then the far field 
intensity can ben expressed by means of the Fourier transform of /, 



If (k x ,k y ) \F{f} (K,k y )\ = 

(5) 

Here, the Jacobi-Anger identity [18] has been used to express the Fourier transform in 
form of Bessel functions J m of the first kind of order m. For a given k y , a diffraction order 
exists, if 3m G Z k y = m2ir /d r . In this case, 



If (k x , k y ) oc 
A comparison with Eq. ()2]) yields 



k x c g j n k x dN oo 

J- . fcfg XI i m Jm(ak x )2ir8(k y -m2ir/d r ) 
, sin 



sin sin 
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(6) 



7/ (A^, m2<K/d r ) = J 2 m (ak x ) x J r (A; x , 0) , (7) 

a formula describing the intensity of the LER grating / given by Eq. (j4]) with respect to 
the undisturbed reference grating r given by Eq. ([1]). The equation states that the impact 
of line roughness does not depend on the roughness period d r , but only on its amplitude a. 
Therefore, the result can be directly generalized to random roughness by Fourier expansion 
of the roughness contours. On the other hand it is sufficient to verify only the sinusoidal 
case, which is much easier to handle, with FEM. In the case of ID random roughness |13j . 
the attenuation factor of Debye- Waller type caused by roughness also was a function only 
of the standard deviation of the random line position or the random line width. In the 
ratio of the diffraction intensities If(k x ,k y ) of the disturbed grating for finite values of k y , 
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to the diffraction intensities l r (k x ,0) of the undisturbed grating, the term l r (k x ,0) cancels 
completely, Eq. ((7j). The ratio only depends on the line roughness. Therefore, this ratio 
can be easily compared to measurements or numerical simulations, even if the values for I r 
do not well represent the measured or simulated values. We present a respective numerical 
comparison in section HI Furthermore, for a given diffraction order n in x, k x = n2ir/d all 
diffraction orders m in y can be summed up: 



^2 If (k x ,m2n/d r ) = l r {k x ,0) ^ J m( ak «>) 

observable m observable m 

< I r (k x , 0) J l Hz) = Ir (k x , 0) . (8) 

ra&L 

Here, "observable m" is a subset of Z referring to the diffraction orders which are ob- 
servable, i.e. reflected back with scattering angles in (—90°, 90°) along the ^/-direction with 
respect to the surface normal. The inequality flB]) states that the undisturbed intensity of 
each order in x is spread through the line roughness over all orders in y-direction for this k x 
value. 

2.B. Line width roughness 

The reflectivity of the absorber structure disturbed by LWR, as outlined in Fig. [TJ^c) is given 
by 

9{X > V) = Pj {X ~ Jd) * ^ ( c- 2a co^/d,) ) 

(9) 

Using the trigonometric angle addition formula and the Jacobi- Anger identity, it follows that 
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c — 2a cos (2-iry/d, 
2 . k x (c — 2a cos (2iry/d r )) 



. ... 2 . k x (c - 2a cos (2ny / d r )) 
exp (-ik x jd) — sin 



exp (— ik x x — ik y y) dxdy 
exp (-ikyy) dy 
(k y ) 
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+ —^J- 2 ^(-1) n J 2n (M vr [5 (k y + 2n2vr/d r ) + 5 (k v - 2n2ixjd r )\ 



2 n =l 

k x c oo 



COS 



+ 2 ^(-l) n J 2n -i (M 7T [5 + (2n - l)2vr/d r ) + 6 (k y - {2n - l)27r/d r )] 

(io; 



Unlike as for LER in the previous section, here it is suitable to discuss the odd and the 
even diffraction orders in y separately to further simplify the expressions. Diffraction orders 
are found at each multiple of the roughness frequency 2n /d r for k y . For the even orders k y , 
3m £ Z k y = 2m(2ir/d r ) , we obtain: 



sin sin ^ 



F{g}(k x ,k y ) (x2vrc— -^^f (-1) H J 2H M . (li; 

sin-*- -f- 

For the odd orders k y , 3m G Z fc y = (2m — l)27r/c? r and ^^0, then 
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cos %p 


sin ^ 


h x c 
2 


sin k *f 


COS %p 


sin %^ 


h x c 
2 



-l) m J 2m _!(M for m>0, (12) 

F{g}{k x ,k y ) oc 2vrc— ^-^^(-1) H+1 J 2H+1 (M for m<0. (13) 

sm * 

This implies that for the even orders 

I g (k x , 2m2n/d r ) = .l\ m (ak x ) x I r (k x , 0) , (14) 

where I g is the intensity distribution of the LWR grating given by Eq. fl9]) an d I r is the 
intensity of the undisturbed reference grating given by Eq. ([1]). For even orders in y, the 
ratio is therefore the same as was found in the case of line edge roughness, cf. Eq. (0). Note 
that for the zeroth order in y, the following approximation is found for small ak x : 

for both LER and LWR. The exponential ratio was also found in the ID stochastic case [13] , 
where the deterministic amplitude a is replaced by the standard deviation of the probability 
density for the line edge. It is worth mentioning that this exponential factor was also found 
independently for small angle X-ray scattering [5]. For the odd orders in y, we find the ratio 

I g (k x , (2m - 1)2tt/4) _ ^m-i ( ak x) 
I g (k x , (2p - l)2n/d r ) ~ 4-i (akx) ' 
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for m,p G N. For LWR, the form factor is no longer the sine function. Therefore, the ratio 
to the undisturbed grating r cannot be discussed as easily as for the even orders. Because of 
the divergence in the odd orders in y, the case k x = must be discussed separately. Setting 
k x = in 

n 

F {/} (K, k y ) = ^ F{exp(-ik x jd) (c- 2a cos (2iry/d r ))}(k y ) (17) 

j=-n 

yields 

F{f}{Q,k y ) =2irN(cS{k y )-a [5 (k y + 2n/d r ) + 5 {k y - 2ir/d r )\) . (18) 

Thus for k x = 0, diffraction orders in the /c^-direction exist only for k y = or k y = ±27r/d r 
with relative amplitudes of c/a, respectively. 

3. Rigorous simulations of light scattering off EUV line masks with sinusoidal 
roughness 

For rigorous numerical simulations of light scattering off EUV masks we use the time- 
harmonic FEM solver JCMsuite [19]. The geometry of a unit cell as schematically depicted 
in Fig. [T] is discretized with a 3D volume mesh. Fig. |2](a) shows a visualization of a mesh for 
typical geometry parameters. Clearly the method also allows for more complex roughness 
models [20]. 

In the simulation, the structure is illuminated with S- and P-polarized plane waves (vac- 
uum wavelength Ao = 13.4 nm) at oblique incidence with an inclination angle of 6°. We 
performed simulations on 3D computational domains with a periodicity in x-direction, 
d = 200 nm, periodicities in y-direction, d r = 200 . . . 1000 nm, and a total thickness of ab- 
sorber stack, Si/Si0 2 capping layer and Mo/Si multilayer coating (60 layers), /i to tai ~ 500 nm. 
The linewidth of the absorber structure was c = 100 nm, with sidewall angle a = 90° and 
top corner rounding radius r roun di ng = 5 nm. We show data for two different roughness peri- 
odicities in y-direction, d r =300 nm and 600 nm, and for six different roughness amplitudes 
a = n x 2.5 nm, where n = 1,...,6. To discretize the vectorial electric field solution to 
the Maxwell light scattering problem we use edge elements of 5 th polynomial order. An 
automatic, rigorous domain-decomposition method is used to separate the essentially ID 
problem of light propagation in the multilayer mirror from the 3D problem in the part of 
the computational domain containing the 3D absorber structure [21]. Figs. |2] (b), (c) show 
visualizations of the EUV light field intensity distribution in two 2D cross sections through 
the 3D computational domain. 

The amplitudes and phases of the reflected diffraction orders in the far field presenta- 
tion are obtained from Fourier-transforms of the 3D near field solution. We note that the 



computational domain is relatively large, with a total volume of 10,000 to 50,000 cubic 
wavelengths. Therefore, the numerical accuracy of the 3D results was checked carefully in 
a previous convergence study and by comparison with a scatterometry experiment [TO} 122]. 
Numerical accuracies with relative errors lower than 0.01% are reached for the central diffrac- 
tion orders with relatively large total power; numerical accuracies with relative errors lower 
than 1% are reached for higher diffraction orders with total intensities down to 10 -6 . The 
results presented here are obtained using a standard multi(8)-processor computer: the typi- 
cal memory usage (RAM) is in the range 10-100 GB, computation times are in the range of 
minutes to tens of minutes. 

4. Comparison simulation vs. calculation 

Figs. |3] (a) and (b) show the diffraction intensities calculated by the FEM simulation. Ac- 
cording to the analytical studies of the previous sections, the major difference between the 
two roughness models should be observed in / (0, ±27r/d r ). This is shown in more detail in 
figs. [3] (c) and (d) where numerical values are presented for the th , 1 st , and 2 nd order in y. 
In the case of LER, a detectable diffraction order for k x = exists only for k y = 0, because 
the Bessel function in Eq. ([6]) vanishes for all other orders in y. In the case of LWR, Eq. ffTB"]) 
states that there are three observable diffraction orders for k x = 0, which are the 0th and the 
±lst order. This is fully consistent with the simulation results shown. In the case of LER, the 
calculated reflectivity in the orders (0, ±1) is below 1 x 10 -6 , whereas there is a detectable 
reflectivity of around 3 x 10~ 3 in the case of LWR, cf. also Figs. |3] (c) and (d). The intensity 
distributions in the even orders as given by Eqs. (J6]) and (flT|) in y are practically identical 
for both LER and LWR as can be seen for the 0th and 2nd order, respectively, in Figs. |3] (c) 
and (d). The analytical results of Eqs. ((7j), (1T4"1) . ( fl6l) are also consistent with the rigorous 
calculation. Fig. H](c) illustrates that, also in the FEM simulation, the ratio of the first order 
roughness induced scatter I (k x ,±27r/d r ) and the undisturbed in-plane diffraction I (k x ,0) 
does not depend on the roughness period d r . The dependence on the roughness amplitude a 
is illustrated in Figs. |4] (a) and (b). For even orders in y, the effect of line roughness is the 
same for both LER and LWR, cf. Fig. [4] (a). Even the ratio for different odd orders as given 
by Eq. (fl6l) has been confirmed by the simulation, cf. Fig. [4] (d). 

Fig. [5] illustrates Eq. ([8]). The undisturbed in-plane diffraction intensities correspond to 
the sum of the disturbed intensities over y. Thus, LER scatters the light in the y-direction 
by conserving the energy for each order in x. 

According to Eq. ([18]) , we obtain 



in the case of LWR. As mentioned before, the Fraunhofer approximation for I r does not 




(19) 
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represent well the real intensities. Particularly, the th order 1(0, 0) is subject to phase effects 
caused by the real 3D structure of the sample as compared to the 2D approximation of the 
Fraunhofer diffraction. These effects are, however, present in measured data and accounted 
for in the FEM calculations. In our case, 1(0, 0) is suppressed in densely patterned fields 
of the photomask as compared to the Fraunhofer approximation. A correction factor can 
be found by normalizing the zeroth order with respect to the first oder in-plane diffraction 
found by the Fraunhofer model and by the FEM simulation as well. 



const 



2/ r (0,0) 
I r (2n/d, 0) + I r (-2ixjd, 0) 



2/ r (0,0) 

I r (2pi/d, 0) + I r (-2pi/d, 0) _ 



2J r (0,0) / fsm(nc/d) 

X 



I r (2-u/d, 0) + I r (-2n/d, 0) \ nc/d 
= 0.468. (20) 

Here, the second term uses the diffraction intensity calculated by the analytical Fraunhofer 
model I r (Eq. fl3])) and I r is the intensity of the undisturbed reference grating calculated by 
FEM. This correction factor can be used to estimate the ratio a/c from the diffraction orders 
for k y = and k y = ±27r/d r at k x = 0: 

°V = const xte^l, (21) 
ci / a (0,0) 

where I g stands for the simulated and potentially measured intensity of the LWR disturbed 
grating. The result for the simulated intensities as used here is shown in Fig. |6j 

5. Conclusion 

We have studied the influence of line edge and line width roughness on the diffraction intensi- 
ties of line gratings, e.g. a scatter field of a photomask. For sinusoidal roughness models, the 
analytical model of Fraunhofer optics leads to a description of the disturbed 2D diffraction 
pattern by means of Bessel functions of the first kind. Also rigorous simulations of the same 
structures were performed using a finite-element-based Maxwell solver. Comparing the ratio 
of the out-of-plane, roughness induced diffraction intensities to the in-plane diffraction of the 
undisturbed structure as calculated by FEM, it could be shown that the ratio fits very well 
the Bessel squared functions, as predicted by Fraunhofer optics. For LER and at least for 
the even off-plane orders caused by LWR, the disturbed diffraction intensities can therefore 
be estimated using the numerical result of an undisturbed line and space grating, and the 
results of the Fraunhofer model without explicitely simulating the line roughness by FEM. 
Some general properties of the 2D diffraction at rough gratings are discussed. In the case 



10 



of LER, the scattering parallel to the lines occurs by conserving the energy in each order 
perpendicular to the lines. For the zeroth order, there are no roughness induced out-of-plane 
scatter intensities for LER and only first order diffraction intensities in the case of LWR. 
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Fig. 1. 2D binary grating without roughness (a). Deterministic 'roughness' 
models, line edge roughness (b) and line width roughness (c). The framed 
boxes indicate the unit cells of the numerical simulation, cf. section [3j 
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Fig. 2. (a) Visualization of parts of a FEM mesh discretizing an EUV line mask 
with sinusoidal LWR, amplitude a = 10 nm and period d r = 600 nm (red: ab- 
sorber, light gray: buffer). Mesh elements discretizing the surrounding vacuum 
and the multilayer mirror below the structure are not shown. Mesh generated 
using the automatic mesh generator JCMgeo. (b) Pseudo-color visualization 
of the electromagnetic field intensity distribution 1(f) on a logarithmic scale 
(colors blue to red correspond to log(J) = — 3 . . . + 1) in a x-z cross section 
through the upper part of the computational domain containing the absorber 
structure, (c) Same, in a x-y cross section through the center of the absorber 
structure. 



14 



-10 -5 5 10 

diffraction order k„ /2nd" 1 



-10 -5 5 10 

diffroction order k, /2nd" 1 



(b) 



(c) 



10"' 

10" 2 

10" 3 

10" 4 

10~ 5 
10" 6 



^ ......... . .......... ^ 


r ® 








r p x' : 




■ \ 


■ f - 




'•A 




if** \ ?\, 


. ici . ™ . i . . . . ! ■ 


ii \ i i H ■ 

I; \./ . , k\ 


-10 -5 ( 


) 5 10 



diffraction order k x /2nd" 1 



10"' 

10" 2 

10" 3 

10- 

10" 5 
10" 6 



(d) 



-x 



,' > 
x >»*\ 



-10 -5 5 10 

diffraction order k, /27Yd~' 



Fig. 3. (a) and (b) Pseudo-color visualization of the result of FEM simulations 
for the structure from Fig. [2] for LER (left) and LWR (right) on a logarithmic 
scale (colors blue to red correspond to log(J) = —5 ... — 1). Numerical values 
are shown in (c) and (d) for the th (red circles), 1 st (blue crosses), and 2 nd 
(green triangles) order in y of Figures (a) and (b), respectively. 
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Fig. 4. Black curves: Bessel squared functions according to Eqs. (JTj), ( JT4l) . ( 1T61) . 
Data points: ratios of the diffraction intensities obtained from FEM simulations 
with respect to the FEM simulation of the undisturbed structure, denoted by 
Iq. LER (red circles) and LWR (blue crosses). The damping ratio for the th 
order in y is shown in (a) for LER and LWR at a = 10 nm and in (b) for LER at 
a = 5 nm both with d r = 600 nm. (c) First order in y for LER at amplitude a = 
10 nm and d r = 300 nm (closed circles) or 600 nm (open circles), respectively, 
(d) Ratio of third and first order in y for LWR at a = 10 nm and d r = 600 nm. 
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Fig. 5. The black polygonal chain connects the diffraction intensities of the 
undisturbed grating (triangles). The red points represent the sum of the 
LER-disturbed intensities over the orders along the y-direction. LER with 
a = 10 nm, d r = 600 nm. 
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Fig. 6. FEM results for k y = ±2ir/d r at k x = in the case of LWR. Shown 
are the ratios of the orders +1 (red circles) and -1 (blue crosses) in y to the 
Qth orc [ er / g (o,0) for 6 amplitudes a at c = 100 nm. The dashed line shows 
the slope a/1 / const obtained from the correction factor derived in Eq. (120"]) to 
account for the difference between J g (0, 0) and 7 g (0,0) as obtained by FEM 
and the Fraunhofer approximation, respectively. 
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